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We present a general mapping between continuous and lattice models of Bose- and Fermi-gases 
in one dimension, interacting via local two-body interactions. For s-wave interacting bosons we 
arrive at the Bose-Hubbard model in the weakly interacting, low density regime. The dual problem 
of p-wave interacting fermions is mapped to the spin- 1/2 XXZ model close to the critical point 
in the highly polarized regime. The mappings are shown to be optimal in the sense that they 
produce the least error possible for a given discretization length. As an application we examine the 
ground state of a interacting Fermi gas in a harmonic trap, calculating numerically real-space and 
momentum-space distributions as well as two-particle correlations. In the analytically known limits 
the convergence of the results of the lattice model to the continuous one is shown. 
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I. INTRODUCTION 

Triggered by the recent successes in the experimen- 
tal realization of strongly interacting atomic quantum 
gases in one spatial dimensional (ID) [U-IH there is an 
increasing interest in the theoretical description of these 
systems beyond the mean field level. Model hamiltoni- 
ans describing homogeneous ID quantum gases with con- 
tact interaction are often integrable by means of Bethc 
Ansatz In practice, however, only a small number 

of quantities can actually be obtained from Bethe Ansatz 
or explicit calculations are restricted to a small number 
of particles. Properties associated with low energy or 
long wavelength excitations can, to very good approx- 
imation, be described by bosonization techniques [Ioj |. 
For more general problems one has to rely on numeri- 
cal techniques such as the density matrix renormalization 
group (DMRG) [Fillip o r the related time evolving block 
decimation (TEBD) [ll GJ]. Both have originally been 
developed for lattice models and thus in order to apply 
them to continuous systems requires a proper mapping 
between the true continuum model and a lattice approx- 
imation. In fact any numerical technique describing a 
continuos system relies on some sort of discretization. 
Here we consider massive bosonic or fermionic particles 
with contact interactions. Only two types of contact 
interaction potentials are allowed for identical, nonrcla- 
tivistic particles, representing either bosons with s-wave 
interactions or fermions with p-wave interactions. Both 
systems are dual and can be mapped onto each other 
by the well-known boson- fermion mapping fl5l . [l6| . A 
proper discretization of ID bosons with s-wave interac- 
tion is straight forward and has been used quite success- 
fully to calculate ground-state [17|, finite temperature 
[ijj, as well as dynamical problems [l9| for trapped ID 
gases. For p-wave interacting fermions a similar, straight 



forward discretization fails however, as can be seen when 
comparing numerical results using such a model with 
those obtained from the bosonic Hamiltonian after the 
boson-fermion mapping. Using a general approach to 
quantum gases in ID with contact interaction (2fj| we 
here derive a proper mapping between continuous model 
and lattice approximation. We show in particular that p- 
wave interacting fermions are mapped to the critical spin 
1/2 XXZ model. By virtue of the boson-fermion map- 
ping the same can be done for s-wave interacting bosons, 
thus maintaining intcgrability in the map between con- 
tinuous and discretized models. As an application we 
calculate the real-space and momentum-space densities 
of the ground state of a p-wave interacting Fermi gas in 
a harmonic trap, as well as local and non-local two par- 
ticle correlations in real space. To prove the validity of 
the discretized fermion model we compare the numerical 
results with those obtained from the dual bosonic model 
as well as with Bethc ansatz solutions when available. 



II. ID QUANTUM GASES WITH GENERAL 
CONTACT INTERACTIONS 

We here consider quantum gases, that are fully de- 
scribed by their two particle Hamiltonian, i.e., the Hamil- 
tonian is a sum of the form 



x o)- 



(1) 



Additionally we require that the true interaction poten- 
tial can be approximated by a local pseudo-potential, i.e. 
it vanishes for Xi 7^ Xj. Since we are in one dimension, 
this leads to the exact integrability of these models in 
the case of translational invariance Q using coordinate 
Bethe ansatz Q . 
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For deriving a discretized Hamiltonian, it is sufficient 
to consider the relative wave function <p(x = x\ — X2) of 
just two particles. The Hamiltonian then reads 



H 



-&l + v(x) 



(2) 



where we have dropped the term corresponding to the 
freely evolving center of mass. 

The continuous two-particle case has been analyzed by 
Cheon and Shigehara 0,1111 • The local pseudo-potential 
V is fully described by a boundary condition on at 
x = 0: Since fulfills the free Schrodinger equation away 
from 0, it must have a discontinuity at the origin as an 
effect of the interaction. Thus we see that 

<j)"{x) x^O 
<5(*)[0'(O+)-0'(O-)] + (3) 
+ S'(x) [0(0+) - 0(0-)] x = 0. 

In the case of distinguishable or spinful [22[ particles both 
singular terms contribute. Due to symmetry, the term 
proportional to the delta function 8 can only be nonzero 
for bosons, while the 8' term exists only for fermions. I.e. 
we have for bosons 



and for fermions 

dt<t>{x) = 



<j)"{x) x^0 
2S{x)<j>'(0) x = 0. 



cj)"(x) x^0 
2<5'(x)0(O) x = 0. 



(4) 



(5) 



In order to get proper eigenstates (i.e. without any sin- 
gular contribution), the pseudo-potential V acting on 
the wave- function must absorb the singular contributions 
from the kinetic energy. Thus the only possible form of 
a local pseudo-potential for bosons is Vb4> = gB8{x)(f)(0), 
while that for fermions reads Vf4> = — gp8' (x)(f>' (0). Note 
that ((/)') is continuous at for bosons (fermions). 
These two possibilities represent the well known cases, 
where the particle interact either by s-wave scattering 
only or by p-wave scattering only, and the interaction 
strength corresponds to the scattering length, respec- 
tively scattering volume, which are the only free param- 
eters left. 

Since all wave functions must have the respective sym- 
metry, we can restrict ourselves in the following to the 
x > sector. We will write 0(0) for \im x ^ + 4>{x) and 
0'(O) for lim a ,_j.o+ 4 ) '{ x )- The above shows that V imposes 
a boundary condition on every proper wave function: 



0'(O) 
0'(O) 



f 0(0) 



bosons, 
fermions. 



(6) 



Eqs.Q and ^ reveal a one-to-one mapping between the 
two cases, i.e., every solution for the bosonic problem 
yields a solution for the fermionic problem with = 
— 4/gp by symmetrizing the wave function and vice versa. 

At this point we emphasize, that boundary conditions 
of the above form are the only ones that are equivalent to 



a local potential [2l|, |23j . While boundary conditions in- 
volving higher order derivatives can be taken into account 
to describe experimental realizations using cold gases in 
quasi ID traps [24j . the necessarily require finite range 
potentials and cannot be described fully by local pseudo- 
potentials. 



III. DISCRETIZATION 

The treatment of continuous gases in one-dimension 
using numerical techniques requires a proper discretiza- 
tion. That is we approximate the two-particle wave func- 



tion 0(x) G L (R) by a complex number 



! (Z), 



where the integer index j describes the discretized rel- 
ative coordinate x = x\ — x<i. We interpret |0^| as the 
probability to find the two particles between (j — h)Ax 
and (j + ^)Ax. In order t o ap ply numerical methods such 
as DMRG or TEBD [H, [Hf efficiently, it is favorable to 
have local or at most nearest neighbor interactions in the 
lattice approximation of the continuous model. It will 
turn out, that the above systems can all be discretized 
using such nearest neighbor interactions only. 

We start with the kinetic term, that can be approxi- 
mated by 



20 7 - + 



y+i 



A.i 



(7) 



In what follows, we will derive two distinct discretiza- 
tions: first for the bosons, where we allow for double 
occupied lattice sites and can therefore use on-site inter- 
actions to reproduce the boundary conditions ©, and 
then for fermions, where double occupation is forbidden 
by the Pauli principle and interactions between neighbors 
are necessary in the lattice model. Note however, that 
both descriptions are equivalent due to the Bose Fermi 
mapping in the continuum limit. 



A. bosonic mapping 

In the lattice approximation the kinetic-energy term, 
Eq.© reads 



dim 



3 > U 

2(0i -0o) • _ n 



(8) 



Thus assuming a local contact interaction only, we find 
for the bosons 



(H4>) 



t/00 



A? 3 > 

Ax 2 3 — u 



(9) 



In order to determine the value of U, we assume, that 
it can be expressed as a series in Ax and evaluate the 
stationary Schrodinger equation {H<j))j — E<j)j = at j = 



3 



0. Rcexpressing <f>\ in terms of </>(0) by means of the 
discrctizcd version of the contact condition (El) 



Ax <f>'(0) ■ 



= -E<f>(0) 



(10) 



we arrive at 



= (H<f>) j=0 - E(j> J=0 
2<t>(0) 



Aa 



Ax^-<j>{0) - ^-Ax z E4>{0) + C(Aa; 3 ) 



- Ecf>(0) - 
1 



A.i 



(11) 



Equating orders gives 



U=^L + O(Ax). 
Ax 



(12) 



The constant term vanishes, since —d 2 (p = E<fi for any 
eigenstate. The higher orders O(Ax) contain E and 
would thus not be independent on the eigenvalue. This is 
perfectly consistent, since discretizations will only work 
a long as the lattice spacing is much smaller than all rele- 
vant (wave) lengths in the system. Thus the lowest order 
in (|12[) is already optimal. There are no higher order 
corrections possible for a general state. 

We can now easily write down the corresponding many 
particle Hamiltonian for the case of indistinguishable 
bosons in absolute coordinates, represented by an inte- 
ger index i and in second quantization: 



H 



E 



-J(a\a i+ i + h.a.) + — a] a J aja . + V.ajai 



(13) 

Here is the bosonic annihilator at site i and V% in- 
troduces an additional external potential in the obvi- 
ous way. So not surprisingly we have arrived at the 
Bose-Hubbard Hamiltonian lattice approximation 

to ID bosons with s-wave interaction. Since Aa; must be 
smaller than all relevant length scales, we are however in 
the low-filling and weak-interaction limits U <C J = ^4-; 
[30| . This does of course not imply that the correspond- 
ing Lieb-Liniger gas is in the weakly interacting regime. 
This result might seem trivial, since we can also directly 
get it by substituting the fiel d operator in the continuous 
model: ty(jAx) i-> 



/Ax 



17|. However, this simple and 



naive discretization does not work in the fermionic case 
we are going to discuss now. 



B. fermionic mapping 

For fermions the kinetic-energy term, Eq.Q reads in lat- 
tice approximation 



dl<t>{x) 



Al? 3 > 1 



Ax 2 





i = i 

3=0 



Due to the anti-symmetry of the wave-function <f>Q must 
vanish, i.e. the simplest way interactions come into the 
lattice model is for nearest neighbors. Thus we write for 
the Hamiltonian 







i>i 
i = i 
j = o 



(15) 



To obtain the value of B we proceed as in the case of 
bosons. As will be seen later on it is most convenient to 
expand B in a series in the following way: 

^ = Ax 2 {b^ + AxB^ +C(Ax 2 )) . (16) 
Now the stationary Schrodinger equation for j : = 1 yields 

= 1 - — Ax- — E + 0(Ax 3 ) + (17) 
3f 2 

B^ + AxB^ + Ax 2 B^ + 0{Ax 2 )^ [l + 0(Ax 3 )] . 

Equating orders results in 

bW = -1, B^ = —, B^ = -E. (18) 
3F 2 

Note that his time the interaction appears only in the 
second lowest order, which can not be described by a 
simple substitution formula. The next higher order con- 
tained in 0(Ax 2 ) does not vanish, but depends again 
on the energy as expected. If we had chosen a straight- 
forward expansion of B instead of (|16[) . the next order 
after the one that introduces the interaction would have 
contained again the interaction parameter: 



B = —■ 



1 



Ax 7 



-J—-± + l +0 {Ax). (19) 
g F Ax gf. 2 



Neglecting this term would therefore introduce a larger 
error than in the chosen expansion (fT6)) . In fact the low 
energy scattering properties would be reproduced only to 
one order less. For the bosons this problem did not occur 
(TT2l . From (fTS)) we read that the optimal result in the 
fermionic case is 



B = — 



1 



1 



Ax 2 \ 1 - 2A-S 



(20) 



The corresponding many-body Hamiltonian for indistin- 
guishable fermions reads 



H = 



\-J(c\ci + \ + h.a.) + Bc\cic\ +1 Ci + \ + V l c\c\ 



(21) 

where now Ci is a fermionic annihilator at site i. Eq. 
(|2"Tj) describes spin polarized lattice fermions with hop- 
ping J and nearest-neighbor interaction B. In contrast to 
the bosonic case, Eg. (TT5)) . where the correct discrctized 
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model could be obtained from the c onti nuum Hamil- 
tonian just by setting ^(x) — > etj/V Ax, wc now sec 
from ([2Tj) and ([20]) that a similar naive and straight- 
forward discretization fails in the case of p-w&ve inter- 
acting fcrmions. 

The failure of a naive discretization of the fermionic 
Hamiltonian becomes transparent if we map this model 
to that of a spin lattice: Using the Jordan- Wigner trans- 
formation 



erf = exp|z7r^ c|q| c i; of = 2c]c 4 - 1 (22) 

Ki 



(f2"Tj) can be mapped to the spin-1/2 XXZ model in an 
external magnetic field 



i 

+ A(of + l)(<rf +1 + 1)) + V;of}, (23) 



where the anisotropy parameter defining the XXZ model 
is A = -l/[l-^l. 

There is an easy way to see that these mappings are 
quite physical by considering the ground states: The re- 
pulsive Bose gas (gu > 0) maps to the repulsive (U > 0) 
Bose-Hubbard model in the super fluid, low filling regime, 
which has an obviously gas like ground state. The 
same is true for the corresponding attractively interact- 
ing (gp < 0) Fermi gas, which maps to the ferromag- 
netic XXZ model which, due to the specific form of the 
interaction parameter in the discretized fermion model, 
Eq. ([2U)) . is always in the critical regime close to the tran- 
sition point (A Aa ~ > ? — 1 + ). A naive discretization would 
have lead to an anisotropy parameter that could cross the 
border to the gapped phase, which is clearly unphysical. 

In the attractive Bose gas, bound states emerge, that 
lead to a collapse of the ground state as it is of course 
also true in the Bose Hubbard model for U < 0. On the 
fermionic side, this collapse can be also observed, as for 
A < — 1 the XXZ model has a ferromagnctically ordered 
ground state, which leads to phase separation in the case 
of fixed magnetization. 

Note that we call the Fermi gas repulsively interacting 
if gp > 0, although B is negative in this case as well, and 
although there exist bound states, who's binding energy 
actually diverges as gp — > + , as is immediately clear 
from the Bose Fermi mapping in the continuous case. 
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FIG. 1: (Color online) Local density distribution of the in- 
teracting Fermi or Bose gas. The (orange) dashed lines show 
results obtained by Bose-Fermi mapping and solving the Bose 
Hubbard lattice model, the (blue) continuous lines correspond 
to the XXZ discretization. The interaction strength qy is 
-51.2, -12.8, -3.2, -0.8, -0.2, -0.05 from the narrow to the 
broad distributions. The solid black lines show the limiting 
cases of free fermions (broad) and infinitely strong interacting 
fermions (narrow, corresponds to free bosons). The calcula- 
tions are done for Ax = One recognizes perfect agreement 
between the fermionic and bosonic discretization approaches. 
Note that both Fermions and Bosons with corresponding in- 
teraction show the same local density, since the quantity is 
invariant under the Bose-Fermi-mapping. 



IV. THE INTERACTING FERMI GAS IN A 
HARMONIC TRAP 

We now apply our method to the interacting Fermi gas 
in a harmonic trap, 

1 N 

j<i 
N i 

+ ( 24 ) 

i=l 

We here chose the trap length to set the length scale. 
For gp — — oo th e sy stem is called a fermionic Tonks- 
Girardeau gas [HI, Hfl [2?J ■ It can be treated analytically, 
since it maps to free bosons under the Bose Fermi map- 
ping. E.g., the momentum distribution is known for arbi- 
trary particle numbers [25j . It is of special experimental 
relevance, since it is equivalent to the density distribu- 
tion measured in a time-of-flight experiment. However 
for intermediate interaction strength numerical calcula- 
tions are required, which we are now able to do. 

First we note, that we now have two options to dis- 
cretized the model. Direct discretization will yield the 
XXZ Hamiltonian, while a Bose Fermi mapping will re- 
sult in the Bose Hubbard Hamiltonian. Both methods of 




k/2jt 



FIG. 2: (Color online) Momentum distribution of the inter- 
acting Fermi gas. Dashed (orange) lines show results via 
the Bose Hubbard discretization, solid (blue) lines corre- 
spond to XXZ discretization. The interaction strength gp 
is -51.2, -12.8, -3.2, -0.8, -0.2, -0.05 from the broad to the 
narrow distributions. Solid (black) lines show the limiting 
cases of free fermions (narrow) and infinitely strong inter- 
acting fermions (broad, calculated from the formula given in 
[25|). The calculations are done for Ax = Again there is 
perfect agreement between bosonic and fermionic discretiza- 
tion. 

course have to produce exactly the same results. 

Fig. [1] shows the spatial density distribution in the 
ground state for N = 25 particles, i.e., 

p(x) = J dx 2 ■ ■ -dx N \<j>{x,x 2 , ■ ■ ■ ,x N )\ , (25) 

which is approximated by the discretized system as the 
diagonal elements of (a\aj). The ground state of the dis- 
cretized system is calculated using a TEBD code and an 
imaginary time evolution, which has already been applied 
successfully to calculate the phase diagram of a disor- 
dered Bose Hubbard model [28[ ■ The interaction strength 
is varied all the way from the free fermion regime to the 
regime of the fermionic Tonks-Girardeau gas. The den- 
sity distribution changes accordingly from the profile of 
the free fermions, showing characteristic Fricdcl oscilla- 
tions, to a narrow Gaussian peak for the fermionic Tonks- 
Girardeau gas. Note that the Bose Fermi mapping does 
not affect the local density, so the curves are the same for 
the corresponding bosonic system. Fe. the density distri- 
bution in the fermionic Tonks-Girardeau regime is iden- 
tical to that of a condensate of non-interacting bosons. 
The curves obtained from the bosonic and fermionic lat- 
tice models are virtually indistinguishable which shows 
that both approaches are consistent. 

The corresponding momentum distribution for the 
fermions, 

Pk(k) = J dk 2 ..- dk N \(f>(k, fc 2 , . . . , k N )\ , (26) 
which is quite different from that of the bosons, is shown 
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FIG. 3: (Color online) single particle density matrix of the 
interacting Fermi gas calculated using XXZ discretization. 
Light regions are positive, dark regions negative. The in- 
teraction strength gp is —51.2,-12.8, and —3.2 (upper row) 
and —0.8, —0.2, and —0.05 (lower row). Remember that the 
cloud size is independent of the particle number towards the 
fermionic Tonks limit (because there is condensation in the 
bosonic picture) while it grow as \/~N for free fermions. 

in Fig. [21 It was obtained from the discretized wave func- 
tion as the diagonal elements of the Fourier transform of 
(ajaj). Again perfect agreement between the bosonic 
and fermionic lattice approximations can be seen. In ac- 
cordance with physical intuition invoking the uncertainty 
relation and Pauli principle, the momentum distribution 
broadens as the real space distribution narrows. While 
for the free particles, real and momentum space descrip- 
tion coincide for the harmonic oscillator the Friedel oscil- 
lations are deformed gradually towards the result for the 
fermionic Tonks-Girardeau gas calculated e.g. by Bender 
et al. [25| . The oscillations that remain in this limit are 
effects from the finite number of particles. They vanish 
as \/N as can be seen from a Taylor expansion in 1/N 
of the expressions given in [25| for the Fermi- Tonks case. 

In Fig. [3] we have plotted the complete single particle 
density matrix 



p(x,y) = 




. . .dx N (j)*(x,X2, ■ ■ ■ )4>{y,x 2 , ■■■) (27) 



for different interaction strength, starting from the 
Fermi- Tonks limit to the case of free fermions. One 
clearly recognizes two small off-diagonal peaks for larger 
interaction strength. The weight of these peaks, which 
are responsible for the oscillations in the momentum dis- 
tribution, Fig. [2] to the remaining part near the diagonal 
is , as can bee seen from analyzing the limiting case nu- 
merically, which can be done for much larger N also. The 
sign of the peaks is positive only if N is odd and nega- 
tive for even N, so the momentum distributions in Fig. 
[3] would show a minimum at k — for all interaction 
strength if N was chosen even instead of 25. 

On first glance it may seem surprising that a map- 
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FIG. 4: (Color online) momentum space distribution of the Fermi gas showing convergence of the method with discretization 
for a) the Fermi Tonks limit, b) qf = —0.8, and c) the free fermionic case. Again in a) and b) dashed (orange) lines show 
results via the Bose Hubbard discretization, solid (blue) lines correspond to XXZ discretization, a) Results are shown for 
Ax = j, | , -jg, yig . As the grid gets finer, both discretization formulas converge to the exact result (black line), b) The 

same discretizations are used as in a) and we again observe convergence of both formulas towards a common limit, which is in 
this case not known analytically. The black lines are those showing up in a) and c) respectively and are for orientation, c) Note 
that in this case there is no sense in distinguishing the two formulas, since implementing U — oo always means excluding double 
occupation of sites by bosons which is immediately equivalent to simulating free fermions. We here only Aa; = j (circles), ~ 
(crosses), (squares) to avoid confusion since the lines converge quite quickly. Although the squares sit perfectly on top of 
the exact result (black lines) they are not spaced densely enough to resolve the Friedel oscillations. This would require a lattice 
that extends across a region in space much larger than N oscillator length where we have chosen to restrict the calculation to 
20 oscillator length to speed it up. 



ping of a continuous, Bethe-Ansatz integrable Hamilto- 
nian such as the Lieb-Liniger model to the non-integrable 
Bose-Hubbard model should produce accurate results. 
However, since the Lieb Liniger gas is dual to p-wave 
interacting fermions, as shown here its lattice approxi- 
mation is equivalent to the spin 1/2 XXZ model, which 
is again Bethe-Ansatz integrable. Furthermore full re- 
covery of the properties of the continuous model can of 
course only be expected in the limit Aa: — > 0. In Fig0]we 
have shown the momentum distribution of p-wave inter- 
acting fermions for decreasing discretization length Aa; 
for three different values of the interaction strength. One 
clearly recognizes convergence of the results as Ax — > 0. 
In the two analytically tractable cases of a free fcrmion 
gas and an the Fermi- Tonks gas the curves approach 
quickly the exact ones. 

As a final application we calculate the real-space two- 
particle correlations in a trap. The corresponding results 
are shown in Fig. [5J Again the (blue) solid lines are 
obtained from the fermionic lattice model and the dashed 
(orange) lines from the dual bosonic model. Due to Pauli 
exclusion g^ifi) = and there is a pronounced dip in 
the g^ near the origin for non interacting or weakly 
attractive fermions, while we see again Friedel oscillations 
for larger inter particle distances. In the dual bosonic 
case the dip is enforced by a strong repulsive interaction. 
As the fermionic attraction is increased, the depth of this 



dip is decreased. There is a smooth transition to the 
perfect Gaussian shape expected for the free bosons in 
the case of strongly interacting fermions. 

Outside the point where the particle positions coin- 
cide both discretization formulas give the same result. 
There is a discontinuity maintaining g^(0) — for the 
fermions, enforced by the symmetry of the wave func- 
tions. It should be noted that this singular jump is not 
reproduced in the dual bosonic model. This is because 
the duality mapping of the discretized models is only 
valid for two particles at different lattice sites and the 
dual bosonic model can only be used to calculate multi- 
particle correlations of fermions at pairwisc different lo- 
cations. 

Finally we note that using the discretization formulas 
(fj"2"j) and (|T5)) one can of course also calculate other many 
body properties like off diagonal order [27} using TEBD 
for larger systems. The method was also used to calcu- 
late out-of equilibrium dynamics for boso nic g ases in the 
repulsive |lj| as well as attractive regime ,29s] . 
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